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We develop a diagrammatic scattering theory for interacting bosons in a three-dimensional, weakly 
disordered potential. We show how coUisional energy transfer between the bosons induces the 
thermalization of the inelastic single-particle current which, after only few collision events, dominates 
over the elastic contribution described by the Gross-Pitaevskii ansatz. 
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Bose-Einstein condensates, initially fascinating by 
themselves, have turned into a playground for a wide 
range of physical phenomena, reaching from condensed 
matter physics to cosmology [HE]- A particularly in- 
teresting subject which recently spurred a lot of experi- 
mental and theoretical activities is the field of coherent 
many-particle quantum transport in disordered potential 
landscapes, due to exceptional experimental control on 
both, the confining potential as well as the inter-particle 
interactions [3j. This allows a detailed scrutiny of the 
hitherto largely unaccessible interplay of many-particle 
[TJ IIHB] and disorder- induced j3j quantum transport phe- 
nomena, at an unprecedented level. In particular, given 
the precise knowledge of the microscopic constituents 
and interactions which define the many-particle eigen- 
states and dynamics, these systems open new perspec- 
tives for an improved understanding of the emergence 
of collective and/or thermodynamic behavior from fun- 
damental dynamical laws, with the thermalization prob- 
lem as one central issue we will here embark on. Al- 
ternatively to various dynamical treatments of thermal- 
ization, in diverse physical contexts [Tllll). we formu- 
late a linear, A^-particle scattering scenario off a disor- 
dered potential. We will see that this lends itself to a 
transparent and physically intuitive understanding of the 
collision-induced thermalization of a quantum system, 
under rather general and experimentally easily accessible 
conditions. Furthermore, in the context of condensate 
dynamics in disordered potentials, our A'^-particle scat- 
tering approach defines a qualitative improvement over 
the wide-spread mean-field (or Gross-Pitaevskii) treat- 
ment, and a first step towards a linear scattering theory 
for interacting many particle systems. 

Specifically, we consider an interacting bosonic gas 
scattering off a three dimensional, weakly disordered po- 
tential. We show that the thermalization is mediated by 
inelastic contributions to the scattering amplitude, that 
amend and rapidly dominate over the strictly elastic, col- 
lective behavior described by the Gross-Pitaevskii equa- 
tion. This establishes a unifying framework for "conden- 
sate depletion" and the formation of a "thermal cloud" , 
as encountered e.g. in [112 , under strictly unitary many- 
particle evolution. The role of the disorder is to ran- 
domize the individual particles' momenta, as necessary 



prerequisite for seeding inelastic collision events. 

While weak particle-particle interactions can be 
treated perturbatively under the assumption that the 
condensate be close to thermal equilibrium [T31 [H] , this 
is not adequate any more in our present situation far 
from equilibrium. We therefore develop a diagrammatic 
theory which involves a non-pcrturbative summation of 
all those contributions which survive the average over the 
weakly disordered potential. 

Let us start with a description of our scattering setup: 
Initially, each atom is prepared in the same single-particle 
momentum eigenstate with wave vector k^, pointing in 
^-direction. Then, the atoms enter a three dimensional 
slab with thickness L along the z-direction, and infinite 
extension in x- and y-direction. Within the slab, the 
atoms experience scattering from a random potential V, 
and collisions due to particle-particle interaction U. On 
exit from the slab, the average spectral flux density J^, 
i.e., the flux of particles with energy E averaged over 
different realizations of the disorder potential ^(r), is 
detected. 

Our microscopic scattering theory starts from an ex- 
pansion of the iV-particle scattering amplitude in powers 
of V and U. Each term in this expansion defines a scat- 
tering diagram and is composed of the following three 
elements: (i) the Green's function Gq for a single par- 
ticle in free space, (ii) scattering of a single particle by 
the disorder potential V{r), and (iii) the two-particle T- 
matrix describing collisions between particles: 
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which, for a short-range interaction potential U and the 
s-wave scattering approximation aski^2 ^ 1, can be ap- 
proximated by |15) : 
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with s-wave scattering length a^, fci — |ki|, fc2 = |k2| 
and fci2 = |ki — k2|. The scattering length together with 
the density po of particles in the initial state defines an- 
other length scale, ^int = l/Snalpo - the mean free path 
between two successive collision events. The probability 
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Figure 1. Example of a ladder diagram describing the propa- 
gation of three interacting particles in a slab with a random 
scattering potential. Pairs of conjugate amplitudes (solid and 
dashed arrows, respectively) undergo the same sequence of 
scattering events (encircled stars) induced by the disorder po- 
tential, at ri, . . . , rg. Due to particle-particle collision events 
(squares), the particles redistribute their energies. Here, solid 
and dashed arrows correspond to disorder averaged single- 
particle Green's functions Q and their complex conjugates, 
respectively. Upon flux detection, one particle is annihilated, 
while the undetected particles are traced over (dots). 



that a third atom is located at distance a^ from the col- 
liding pair is assumed to be small, a^po ^ 1, such that 
three-particle collisions can be neglected jlB| . 

The A^-particle scattering amplitude obtained by all 
combinations of the above building blocks (i - iii) defines 
the final state |/) = S\i), where \i) = {Nki) is the ini- 
tial state with all N particles in mode k^. The measured 
flux Je — {f\JE\f) is derived from the single-particle 
observable Je which annihilates one particle, with the 
remaining A^ — 1 particles traced over. To end up with 
a statistically robust quantity, one finally needs to aver- 
age over different realizations of the disorder, which we 
assume Gaussian distributed, with mean V(r) — 0, and 
correlation 
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where ^dis is the mean free path for scattering off the dis- 
order potential. Under the further assumption of weak 
disorder, fc^dis S> 1, only so-called ladder diagrams |17j . 
where the amplitude and the conjugate amplitude un- 
dergo the same sequence of scattering events, survive the 
disorder average. Go is replaced by the average single- 
particle Green's function 
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where Jce = vE + «/2£dis, in units of h^ /2m. 

A diagram contributing to the average flux is then 
constructed as follows: we take one diagram for the A^- 
particle scattering amplitude, another one for the con- 
jugate amplitude, group them together into a ladder di- 
agram, detect one of the outgoing particles, and trace 
over the other ones. A typical example is shown in Fig. IT] 
Among the A^-particle ladder diagrams thus constructed. 



we neglect all those where two particles which interacted 
once meet again. Alike the neglect of non-ladder dia- 
grams, this approximation is valid for fc^dis ^ 1, and al- 
lows us to trace over the undetected particles after their 
interaction with the detected particle, as shown in Fig.jT] 
Finally, we assume that at least one disorder scattering 
event occurs between two collision events. This is justi- 
fied if £int ^ ^dis, which, for a realistic scenario, is easily 
accessible by adjusting £disj i-e., the disorder strength, 
accordingly. 

Under these assumptions, any diagram contributing to 
the particle fiux is composed of three building blocks, 
see Fig. [2J The first one. Fig. l2Fa), denotes scattering 
of a single particle off the disorder potential at ri, and 
subsequent propagation to the next scattering event at 
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where Vyi — '^\ — '^i- In the second building block. 
Fig. [2Fb), one pair of amplitudes (solid lines) exhibits 
a particle-particle collision event, whereas the other pair 
(dashed lines) does not experience a collision. Conse- 
quently, the energies E\ and E2 of both particles remain 
conserved (otherwise the solid and dashed amplitudes 
could not be grouped together). The diagram Fig. [2jb) 
hence represents a nonlinear elastic scattering contribu- 
tion, which we denote by gEi;E2 (ri, r2 7 ra). Furthermore, 
one can show that, if one neglects the second order term 
ki2al in the two-particle scattering amplitude, Eq. ([2|, 
this diagram is equivalent to a diagram obtained from the 
stationary Gross-Pitaevski equation [18 . In contrast, the 
diagram Fig. [2Fc) represents inelastic scattering events, 
not accounted for by the Gross-Pitaevski equation: It 
describes a collision between two particles, where the en- 
ergies of both particles change from Ei and E2 to E^ 
and E4 = E1+E2 — E3, respectively. The weight of such 
processes is given by: 
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Note that Eq. (|6| is quadratic in the small parameter 
fls, i.e., inserting the first order contribution to i(fci2) in 
Eq. (I2]) suffices. The first factor 2 in Eq. ^ originates 
from the fact that the solid and dashed incoming ampli- 
tudes can be grouped together in two different ways. It 
can be shown that this accounts for fluctuations of the 
atomic density inside the disordered slab |19^. The sum 
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Figure 2. The three building blocks from which all ladder di- 
agrams (see Fig.[l| are constructed, (a) Single-particle prop- 
agation in the disorder potential, see Eq. ([5|. (b) Nonlinear 
elastic scattering gEi;E2- (c) Nonlinear inelastic scattering 
/bi,E2;B3i see Eq. ^. 



over k represents the trace over the undetected particle. 
The factor 1/2 in front of the second sum indicates that 
this sum is taken over the subspace of symmetrized two- 
particle states |ki,k2). 

We can now write down a nonlinear integral equation 
for the average particle density at energy E: 

lE{r) = Io{v)S{E-E,) + Jdv'PE{r',v)lE{r')+ (7) 
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which upon iteration produces all possible combinations 
of the above three building blocks, non-perturbatively 
in the collision contributions Sei-.e and fEi,E2\E- Here, 
/o(r) = /Ooexp(— z/^dis) denotes the density po of par- 
ticles in the initial mode, attenuated by the propaga- 
tion to position z inside the slab. Furthermore, we 
have employed a contact approximation for the collision 
terms, i.e., gBi;B(ri,r2,r) ~ 5{Yi~Y)5{Y2-r)gE^-E, with 
gEv,E = /dridra gEi-EijCi,T2,Y), and similarly for /. 
This is justified since we assume ^mt ^ ^dis, and hence 
the spatial transport of particles between two points r' 
and r is dominated by the single-particle propagator 
Pe{'^' , r). With this approximation, Eq. (l6| implies, with 
fc^dis > 1: 
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for E2 > El. The expression for QEi-.e can be calcu- 
lated in a similar way from diagram Fig. I2|b). Alterna- 
tively, it can be extracted from Eq. (Is]) and one of the 

two conditions V^5Bi:B2 = ^ /o ' ^ ^-^ VEfEi,E2;E 

or {El + E2)VE^gE,;E2 = - Jo^'^'^" dE 2Ey/EfE,,E2;E. 
which guarantee conservation of the particle and energy 
flux, respectively, in Eq. ([7]). In other words, inelastic 
scattering goes along with a corresponding reduction of 
the nonlinear elastic component. To this end, it is crucial 
to keep the second order term in Eq. ([2| in the expression 
for g, since the first-order term (i.e. the result predicted 



by the Gross-Pitaevski equation) vanishes within the lad- 
der approximation |18j . 

With these premises, we can now infer the total av- 
erage flux J{z) ~ J dEJE{z) at position z E [0; i] 
within the slab, with the energy-dependent flux Je{z) — 
\EIe{z), in units of the incident flux Jo = \/EiPo. Here, 
Ie{z) is obtained numerically via iterative solution of 
Eq. Q. Fig. Isla) shows the result for a slab thickness 
h — £/^dis = 50, and weak interaction €dis/Ant = 1/250. 
J{z) exhibits the characteristic linear decay of diffusive 
(or Ohmic) transport [20", and equals the linear flux 
which is obtained from Eq. ([7]) when setting Os = 0. 
This is due to the condition £i,it S> -^dis, and the cor- 
responding contact approximation mentioned above, to- 
gether with the fact that, for the 3D white-noise po- 
tential ([3]), £dis is independent of the particle's energy 
|20) . In contrast to the linear case, however, J{z) splits 
into an elastic and an inelastic component, defined by 
Je{z) = jf\z)5{E - E,) + J^e'''^\z). We see that, in 
spite of the weakness of the interaction, the inelastic com- 
ponent rapidly dominates as the particles penetrate into 
the slab. This can be explained by the large number b^ 
of (disorder) scattering events required to traverse a slab 
with thickness b. The expected number of two-body col- 
lision events correspondingly scales as fe^^dis/Ant- By the 
same argument, three-body collisions can be neglected if 
fljPo ^ Ant/(^dis&^)- Note that the continuous emergence 
of an inelastic component of the fiux, as described by our 
present, microscopic and strictly unitary treatment, is 
tantamount to the formation of what is colloquially called 
a "non-condensed fraction" or "thermal cloud" , since an 
A'^-fold product of a single-particle state (as required from 
the formal definition of a condensate via the stationary 
one-particle density matrix [21 ) with fixed total energy 
implies fixed energies also for the individual particles. 

The normalized energy distribution J]j"° (z) of the 
inelastic component is shown in Fig. IsFb). We see 
that, starting out from a distribution which is cen- 
tered around Ei after the first inelastic event, the en- 
ergy distribution approaches a Maxwell-Boltzmann dis- 
tribution j|f^^ = 4£;exp(-2£:/^i)/^f, with the av- 
erage energy (or "temperature") fixed by the incident 
particle energy E-^. Note that due to the strongly re- 
duced density, collisions become very unlikely towards 
the end of the slab, and hence the spectral inelastic flux 
densities between slab positions z = L/3 and z — L 
are only slightly altered. Fig. IsFb). This also mani- 
fests in Fig. IsFa), by the saturation of the difference 
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||(j(incl) _ j(MB))(^)|| ^ JEiJdE (j^""^'^ - J),"^"^j (z) 

between the inelastic and Maxwell-Boltzmann flux den- 
sity, as a function of position in the slab. 

In summary, we formulated a microscopic transport 
theory for interacting bosons propagating in a ran- 
dom potential. We showed that the disorder-averaged 
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Figure 3. (Color online) (a) Flux density J vs. position z 
in a slab of thickness b — L/^dis = 50, for ^dis/^int = 1/250. 
Elastic (green, dotted) and inelastic (blue, dash-dotted) com- 
ponent add up to the total, interacting many particle flux 
density (black, solid), which coincides with the linear flux 
(red, dashed). The black dash-double-dotted line indicates 
the deviation of J{z) from a Maxwell-Boltzmann flux den- 
sity (y-scale on the right), (b) Normalized energy distribution 
j]j"° (z) of inelastically scattered atoms at z = (blue, dash- 
dotted), z = L/3 (green, dotted) and z — L (red, dashed) 
within the slab; otherwise the same parameters as in (a). 
The distribution fEi,Ei;E after one single inelastic scattering 
event (gray-shaded area), Eq. ([8|, evolves into a Maxwell- 
Boltzmann distribution (black, solid) with average energy Ei, 
upon penetration into the slab. 



tionary energy distribution is known to be given by J^ 
[23j . We note that rigorous derivations of a nonlinear 
quantum Boltzmann equation similar to the coUisional 
terms in Eq. ([t]) have been attempted recently [211 Hi] , 
though in the absence of a random potential. It is pre- 
cisely the presence of the latter, however, which allows 
for a rigorous quantification of the regime of validity of 
Eq. (It]), in terms of the parameters fc, ^dis and ^int, in our 
present treatment. 

Let us conclude with a discussion of possible exten- 
sions of our theory: Since the main idea - neglect of all 
but ladder diagrams for weak disorder - is not restricted 
to stationary scattering processes, we expect that the 
present theory can be extended to time-dependent sce- 
narios such as, e.g., expansion of an initially confined 
condensate in a random potential |26j . Furthermore, re- 
laxing the contact approximation (assuming l-mt 3> ^dis) 
for the collision terms allows to enter a regime of stronger 
interactions where, e.g., repulsion or attraction between 
particles will affect the spatial density profile. Finally, 
by means of crossed diagrams, one can study interfer- 
ence phenomena like coherent backscattering |i27ii, and 
clarify how these are modified by interactions. 
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single-particle density matrix relaxes to a stationary 
state which, after only few (6^-^dis/^int ~ 10) collision 
events inside the scattering region, coincides with a 
thermal Maxwell-Boltzmann distribution with "temper- 
ature" given by the incident particles' energy. 

For interacting particles with confinement rather than 
disorder, an analogous result was derived from a random 
matrix argument, under the constraint 1/k < a^ ^ pg 
[7]. Observe that the assumptions required for our the- 
ory, ke^is > 1 (weak disorder), ^int = l/(87ra^po) > ^dis 
(collisions less frequent than disorder scattering), and 
the s-wave scattering approximation, can be summarized 
as Ug ^ l/k <C ^dis "C ^int- Since k and ^dis may 
be adjusted by choosing the initial energy and the dis- 
order potential appropriately, our results hold for typi- 
cal experimental parameters of ultracold bosonic gases, 
(10^^ - 102^)m-3 and a, w (10^* - 10-^)m [22], 



corresponding to ii^t ~ (4 x 10~^ — 4 x 10~'^)m. We ver- 
ified, for a wide range of optical thicknesses b and ratios 
^dis/-^int, that Eq. ([7| indeed provides unique stationary 
solutions and predicts full quantum thermalization of the 
flux density if the number of inelastic collisions exceeds 
6^£dis/Ant ~ 10. In fact, the collision terms in the trans- 
port equation Hj exactly reproduce Boltzmann's kinetic 
equation for a gas of classical particles, for which the sta- 
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